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Abstract 

An  anti-noise  problem  on  a  finite  time  interval  is  solved  by  minimization  of  a  qua¬ 
dratic  functional  on  the  Hilbert  space  of  square  integrable  controls.  To  this  end,  the  one- 
dimensional  wave  equation  with  point  sources  and  pointwise  reflecting  boundary  conditions 
is  decomposed  into  a  system  for  the  two  propagating  components  of  waves.  Wellposedness  of 
this  system  is  proved  for  a  class  of  data  that  includes  piecewise  linear  initial  conditions  and 
piecewise  constant  forcing  functions.  It  is  shown  that  for  such  data  the  optimal  piecewise 
constant  control  is  the  solution  of  a  sparse  linear  system.  Methods  for  its  computational 
treatment  are  presented  as  well  as  examples  of  their  applicability.  The  convergence  of  dis¬ 
crete  approximations  to  the  general  optimization  problem  is  demonstrated  by  finite  element 
methods. 
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1.  Introduction 

This  article  s  optimization  problem  is  motivated  by  a  variant  of  active  noise  control  in 

acoustics.  In  physical  terms,  we  consider  a  wave  guide  of  finite  length  with  the  following 
properties: 

(i)  the  evolution  of  sound  in  the  guide  is  appropriately  described  by  plane  waves  prop¬ 
agating  along  its  axis  without  dissipation 

(ii)  the  sound  in  the  guide  originates  from  initial  conditions  and  sources  with  small 
support  located  within  the  wave  guide  or  at  its  front  ends 

(iii)  the  system’s  total  behaviour  results  from  superposition  of  mutually  independent 
sound  waves  that  are  generated  by  the  sources  and  the  initial  conditions 

(iv)  at  the  front  ends  of  the  guide  incident  waves  are  partially  absorbed  and  partially 
reflected  without  alteration  of  shape. 

The  initial  conditions  at  time  t  =  0  are  given  and  the  behaviour  of  some  of  the  sources 
(offending  but  uncontrollable  sources)  is  known.  Our  control  problem  is  to  determine  the  be¬ 
haviour  of  the  other  sources  (controllable  sources)  such  as  to  reduce  the  noise  by  destructive 
interference  in  prespecified  regions  of  the  wave  guide. 

Most  of  the  above  assumptions  are  used  in  the  literature  on  anti-noise  problems  in  ducts 
(e.g.  [Sw],  [TN],  [FW],  see  also  [MI]  p.467).  Furthermore,  in  acoustics  it  is  common  practice 
to  consider  harmonic  sound  fields  with  time  dependences  of  the  form  exp(— iw<).  In  contrast 
to  this,  the  present  paper  works  with  non-harmonic  sources  and  fields;  for  the  time  domain 
model  considered  here,  a  convenient  choice  for  the  states  and  sources  are  elements  in  Hilbert 
spaces  of  real-valued  functions. 

The  basic  time  domain  model  in  linear  acoustics  is  the  wave  equation  for  the  velocity 
potential  <f> 

(1.1.1)  </»te(*,a^)  -  c^<^xz(fja:)  = /(f,x) -f  u(f,r),  0<t<T,  0<x<L 

(1.1.2)  <f>(0,x)  =  (f>(i{x),  <^t(0,x)  =  ^o(aj),  0<x<L 

with  initial  conditions  ^OjV’o  :  (OjT)  — ♦  R,  the  constant  c  >  0  being  the  speed  of  sound.  As 
is  customary  in  anti-vibration  models  (e.g.  [Sw],[NCEBl],[NCEB2],[P])  we  consider  point 
sources  that  are  represented  by  6  distributions 

n  m 

(1.1.3)  ='Y^fi{t)6^,{x),  u(f,x)  =  ^u,(f)6,,.(x),  0  <  f  <  T, 

1=2  i=2 

at  distinct  points  r/j  G  (0,L)  with  time  dependencies  G  L^{0,T). 

According  to  (iv),  the  specific  acoustic  impedances  of  the  boundary  surfaces  are  real 
constants  Co,Ci  ^  0-  The  interaction  of  a  surface  zmd  the  sound  is  such  that  the  quotient 
p/t’«n  of  the  sound  pressure  p  imd  the  incident  velocity  Vin  at  the  surface  is  equal  to  pcQ 
([MI]  p.259f,  p  denotes  the  equilibrium  density  within  the  duct).  In  terms  of  the  velocity 
potentiad,  since  p(t,z)  =  p(f>t{t,x),  v{t,x)  —  — <^x(f,x),  this  leads  to  the  boundary  conditions 

<^t(<,0) -cCo<^x(t,0)  =  0, 

(^t(<,  L) -f- cCi  <^x{t)  T)  =  0, 
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We  include  sources  at  the  boundaries,  assuming  that  sources  and  reflection  combine  by 
superposition.  For  exan'  the  model’s  boundary  at  a;  =  i  could  be  an  artificial,  non¬ 
reflecting  one  (Cl  =  1)  ’  udependent  waves  entering  from  x  >  L.  Let  us  consider  the 

case  of  an  offending  soui .  c  at  i  =  L  and  a  controllable  source  at  x  =  0, 

<;i<(<,0)  -  cCo«ji*(<,0)  =  ^  ^ 

(1.1.4)  2c  ^>0 

+  cCMt.L)  = 

/i,ui  €  X^(0, T).  The  factors  (1  Ci)/2c  have  the  effect  that  the  waves  generated  at  the 
boundaries  and  in  the  interior  have  equal  amplitudes  relative  to  fi{t),Ui{t). 

p{c^<j>x{t,x)‘^  -I-  <i)t{t^x)^)l2c^  being  the  acoustic  energy  density  at  time  i  ([MI],  (6.2.15)), 
we  define  the  performance  index 

(1.3)  J(iti,...,Um)  =  /(^  j  q{t,x){c^<t>x{t,xf  +  4>t{t,xf)dx +  Y^ri{t)ui{tf)dt 

0  0 

with  design  parameters  ri{t)  >  0,  i  =  l,...,rn.  The  weighting  function  q{t^x)  >  0  is 
non-zero  on  those  subintervals  of  (0,  L)  where  the  acoustic  energy  is  to  be  reduced. 

The  optimal  control  problem  in  terms  of  (1.1)  is:  Given  initieil  conditions  (f>oii’o  and 
functions  fi  €  L^(0,  T),  i  ~  1, . . . ,  n,  find  functions  Ui  €  L^(0,  T),  i  =  1, . . . ,  m  that  minimize 
the  functional  J  wherein  <f)  is  the  solution  to  (1.1).  Inserting  the  solution  operators  of  the 
system  into  the  cost  functional  J  gives  a  quadratic  functional  in(ui,...,Um)  whose  unique 
minimizer  abstractly  is  represented  by  use  of  the  adjoint  of  the  operator  that  maps  the  control 
functions  to  the  waves  they  generate.  The  purpose  of  the  present  work  is  the  construction 
of  approximate  optimization  problems  that  are  amenable  to  numerical  computations. 

To  this  end,  equation  (1.1)  is  substituted  by  a  first  order  system  according  to  the  fac¬ 
torisation  dtt  —  c^dxx  =  {dt  —  cdx)(dt  4-  cdx)-  The  components  of  that  system  are  the  two 
components  of  the  velocity  potential  <}>  that  are  travelling  in  opposite  directions.  The  delta 
inhomogeneities,  generally  speaking,  involve  non-classical  but  wellknown  concepts  of  solu¬ 
tion  (method  of  transposition,  e.g.  [L]).  Yet,  for  the  boundary  conditions  used  in  this  article, 
no  directly  quotable  reference  that  would  cover  the  wellposedness  of  the  models  in  section  2 
was  found.  Instead  of  following  the  general  theory  of  [LM],  for  the  one- dimensional  systems 
at  hemd,  the  existence,  regularity  and  uniqueness  of  weeik  solutions  (in  a  sense  analoguous 
to  [L])  can  be  proved  directly  by  application  of  classiceil  methods  (d’Alembert’s  formula, 
Duhamel’s  principle,  energy  estimate)  in  the  sense  of  distributions. 

The  regularity  of  the  initiad  conditions  and  somces  in  section  2  is  taylored  precisely  to  the 
discretized  minimization  problem  of  section  3.  By  restriction  to  sources  that  are  piecewise 
constant  (in  time)  and  to  initial  conditions  that  eire  first  order  splines  (in  space),  the  problem 
cam  be  formulated  in  terms  of  matrices.  This  is  due  to  the  fact  that  -  according  to  the  simple 
process  of  propagation  and  reflection  -  the  value  of  the  cost  functional  is  determined  by  the 
two  components  of  (f>  at  the  points  Xj  =  jcAt  at  discrete  times  ik  =  kAt.  Thus,  for  this  type 
of  sources  amd  initial  conditions,  the  problem  with  finite  time  horizon  (T  <  oo)  is  reduced 
to  the  solution  of  a  finite  (possibly  large)  system  of  linear  equations. 
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The  coefficient  matrix  of  this  system  is  sparse.  We  present  numericaJ  methods  that  avoid 
the  zero  coefficients.  To  demonstrate  their  computational  fezisability  and  to  illustrate  some 
characteristic  details  of  the  solutions,  four  examples  are  given. 

Finally,  twice  differentiable  initial  conditions  are  interpolated  by  linear  splines  zuid  square 
integrable  sources  are  projected  onto  subspaces  of  piecewise  consteuit  functions.  It  is  shown 
that  the  corresponding  optimal  controls  converge  to  the  solution  of  the  original  problem. 
The  proof  rests  on  finite  element  methods  for  elliptic  problems. 

As  is  generally  the  case  for  wave  equations  (cf.[DLP],[DW]),  the  performemce  of  the 
controls  is  sensitive  to  time  delays.  This  feature  affects  the  open  loop  solution  of  the  present 
paper  as  well  as  the  feedback  solution  as  in  [BKSWj.  It  should  be  noted  that  both  approaches 
lead  to  acausal  controls  in  the  sense  that  the  behaviour  of  the  offending  sources  at  times  later 
than  t  is  used  to  determine  the  optimal  control  action  at  time  t.  This  is  necessary  because  of 
the  global  nature  of  the  performance  index  and  the  unrestrictive  choice  of  admissible  control 
functions.  In  the  exjimples  we  can  observe  the  correlation  of  acausal  actions  with  the  finite 
speed  of  propagation. 

For  typical  boundary  materials  that  are  of  practical  interest,  assumption  (iv)  is  not  valid 
([BG],[BPS]).  However,  modelling  and  approximating  the  smearing  of  incident  pulse  waves 
by  the  boundaries  (see  [MI]  p.264f)  would  destroy  the  sparsity  in  the  matrix  of  the  discretized 
system  and  thereby  drastically  increase  the  computational  requirements. 

On  the  other  hand,  we  believe  that  explicitly  computable  solutions  to  minimization  prob¬ 
lems  for  pointwise  reflecting  boundaries  as  in  (iv)  provide  detailed  insights  into  the  mecha¬ 
nisms  of  the  optimal  control  of  waves  in  enclosures. 

Notation.  L'^{a,h\X)  is  the  Lebesgue  space  of  square  integrable  functions  on  the  in¬ 
terval  (a,  6)  with  values  in  X.  For  real-valued  functions,  A"  =  R,  we  simply  write  L^{a,b). 
C([a,  6];X)  denotes  the  space  of  X-valued  functions  that  are  continuous  on  the  closed  in¬ 
terval  («,  6).  H^{a,b)  is  the  Sobolev  space  of  absolutely  continuous  functions  h  :  (a,  6)  — »  R 
with  distributional  derivative  h'  G  L^(a,b),  normed  by  ||/i||i  =  (||/i|p  +  (||  •  ||  being 

the  JD^-norm).  As  usual,  X^  =  {(xj , . . . ,  x*)|x,-  €  X,  i  =  1, . . . ,  fc}.  In  a  Hilbert  space  X  the 
inner  product  will  be  denoted  by  (’,’)x-  We  use  the  abbreviations  h(^'*')  and  h(^~)  for  the 
right  and  left  hand  limit  of  the  function  h  at  the  point  ^  €  R.  For  am  interval  /  C  R  we 
write  x-f  for  the  characteristic  function  of  J,  i.e.  xA^)  =  1  for  t  G  /,  X-^(0  =  0  for  <  G  R  \  /. 


2.  The  first  order  system 

Any  unidirectional  linear  wave  is  the  sum  of  two  components  propagating  in  positive  and 
negative  direction.  The  insertion  of  the  sum  of  a  function  of  (t  —  x/c)  and  a  functi  m  of 
(<  -f-  x/c)  into  (1.2)  shows,  that  the  boundary  conditions  (1.2)  couple  the  two  compc.'^ents 
by  pointwise  reflection  with  reflection  coeflBcients 


0,1. 
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0  <  X  <  L,  0  <  i  <  T; 


Therefore  we  consider  the  system 

+  c<f>2'(i,x)  =  0, 

=  0» 

<f>-^(t,0)  =  Jio<^-(t,0),  <^-(t,L)  =  Rj</>-^(t,L),  0<i<T; 

<;f>'^(0,x)  =  </>^(x),  <^~(0yx)  =  <^q(x),  0  <  X  <  L. 

This  first  order  system  is  a  generalization  of  the  homogeneous  second  order  equation 

(2.4)  <j>tt{t.,x)  —  (?<i>xx{t^x)  =  0,  0<x<L,  0<t<T 

with  initial  and  boundary  conditions  (1.1.2),  (1.2)  in  the  following  sense:  If  are 

functions  satisfying  (2.1)  -  (2.3)  with  initial  conditions 

I  I 

0  0 

then  <f>  =  (f)'^  +  (l>~  IS  3i  solution  of  (2.4),  (1.1.2),  (1.2)  in  the  sense  of  distributions  [S2];  this 
can  be  seen  by  the  application  of  the  commuting  product  {dt  +  cdx)(dt  —  cdx)  to  <f>'^  +  (f>  , 
and  by  verification  of  the  initial  and  boundary  conditions.  Classical  solutions  that  are  twice 
continuously  differentiable  evolve  only  if  the  initial  conditions  satisfy  certain  regularity  and 
compatibility  requirements. 

2.1.  Wellposedness.  In  context  of  the  present  anti-noise  problem  it  is  adequate  to  work 
with  continuous  initial  data,  because  typical  initial  conditions  are  either  silence,  or  sound 
that  has  been  generated  by  sources  that  axe  driven  by  L^{0,T)  functions.  Accordingly,  we 
define  the  space  of  initial  conditions 

y  =  {(<^+,<^-)  €  H\0,Lf\<f>-^{0-^)  =  Ro<f>-{0+l(f>-{L-)  =  i?i<^+(L-)}, 
endowed  with  the  if^(0,L)^-norm. 

Definition  2.1.  Given  (<^J,  <^^)  G  V,  a  pair  <f>~)  of  functions  :  [0,  T]  x  (0,  L) 

R  is  called  a  solution  of  system  (2.1)  -  (2.3)  if 

(i)  <f>'^i0,x)  =  <f>i{x),  <f>~{0,x)  =  <f)Q{x),  0  <x  <  L; 

(ii)  the  mapping  [0,  T]  —>■  V,  defined  by  <  i-»  ((^"’’(f,  •)>  <^~(^»  ■))>  C([0,T];  V); 

(hi)  for  all  X  E  (0,L),  the  functions  [0,  T]  — >  R  defined  by  <  t-v  <^+(<,x),  t  v->  <f>~(t,x)  are 
absolutely  continuous  on  [0,T]; 

(iv)  for  all  i  €  (0,  T),  (2.1 )  holds  for  almost  all  x  e  (0,L). 

Theorem  2.1.  Given  (<f>o,<^o)  ^  system  (2.1)  -  (2.3)  has  a  unique  solution.  With 
£(t)  =  max{£  G  <  ct/Lj  the  solution  is  given  by 

for  £(i)  =  I  odd: 

x)  =  Ro{R\ Ro)^ <f>Q  {—{£  -  l)L  —  X  4-  cf)  -f-  (iZo-Ri)”^  +  1)-^  +  x  —  cf), 
(f>~(t,x)  =  R\{RoR\)~ +  1)-^  —  X  —  c<)  -f  (i2iilo)'^<^o^(— (^  +  1)L  -|-  x  -f-  ct)\ 


(2.1) 

(2.2) 

(2.3) 
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/or  i{i)  —  i  even: 

x)  =  (RoRi)^  <^i(eL  +  x-ct)  +  RoiRiRo)Uo(-^L  -x  +  ct), 
x)  =  (RiRo)i(^- {-ei  +  x  +  ct)  +  i2i(i2oi2i )^o'((^  +  2)L-x-  ct). 


These  formulae  are  understood  according  to  the  convention  <}>q{t)  =  =  0  if  r  € 

^  \  (0,-C'),  so  that,  for  every  {t,x)  G  [0,  T]  x  (0,  L),  at  most  one  of  the  two  terms  in  each 
formula  is  non-zero.  At  the  points  x  =  —£L  -f  ct  and  i  =  (^  -}-  1)L  —  ct,  where  the  values  of 
and  (f)~  are  left  open  by  this  convention,  these  values  are  determined  as  the  left  or  right 
hand  limits,  which  coincide  due  to  the  boundary  conditions  in  V. 

Proof.  As  to  uniqueness,  the  classical  form  of  an  energy  estimate  (eg.[S2]  p.299)  is  not 
applicable  here,  because  the  initial  conditions  and  solutions  in  Def.  2,1  axe  not  continuously 
differentiable.  However,  for  any  solution  of  (2.1)  -  (2.3)  consider  the  function 

L 

E(t)  =  J  {<l>'^{t,x)^  +  <f>~{t,x)^)dx. 

0 

Because  of  (ii),  (iii),  E  can  be  viewed  as  a  distribution  whose  derivative  is  given  by 


L 

E'{t)  =  j {<f>'^{i,x)^  -f  <f>~(t,xf)tdx 

0 

(see  [Z],  2.8).  TVansforming  the  integrand  using  (iv)  and  integrating  x)^  <f>~{t,  a:)^)* 

from  0  to  T  we  get,  (f>'^{t,  •y^,(f>~{t,  •)^  being  absolutely  continuous  (eg.  [HS],  (18.16)), 

E\t)  =  c[<^+(<,0+)2  -  <^+(t,T-)2  -h  (l>-{t,L-y  -  <f>-{t,0+)^] 

=  c[{Rl  -  l)<l>-{U0+)^+{R]  -  <  0. 


This  shows  that  the  derivative  of  is  a  non-positive  function;  therefore  ([Si], Chap. IV)  E{t) 
is  a  decreasing  function.  Thus,  given  {(i>t^<f>o)  €  V,  the  difference  of  two  solutions  of  (2.1)- 
(2.2)  is  a  solution  with  zero  initial  conditions  i.e.  E{0)  =  0.  Consequently,  for  the  difference 
of  two  solutions,  E{t)  =  0,  <  G  [0,  Tj.  This  implies  the  uniqueness  of  the  solution. 

The  formulae  for  x),  x)  describe  the  movement  of  the  two  components  due  to 
propagation  and  reflection  for  the  (£  -I-  l)th  cycle  of  complete  reversion  during  iLjc  <  t  < 
{i-{-l)L/c  (cf.  section  2.2).  In  light  of  the  convention  stated  above,  it  is  quite  straightforward 
to  verify  that  the  given  pair  satisfies  (i)  -  (iv)  of  Definition  2.1,  □ 


With 


H  :=  L'^{0,T;L\0,Ly) 


we  have  the  following  immediate  consequence  of  the  formulae  in  Theorem  1.1. 
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Corollary  2.1.  The  operator  S  :  V  H  that  maps  any  initial  condition  in  V  to  the  spatial 
derivatives  of  the  corresponding  solution  of  (2.1)  -  (2.3),  =  ((;/»+(<,  ),  •)), 

where  ,<f>~  are  given  in  Theorem  2.1,  is  a  bounded  linear  operator. 

With  regard  to  sources  at  the  boundaries  we  consider  the  homogeneous  first  order  system 
(2.1)  with  zero  initial  conditions 

(2.5)  ^+(0,x)  =  0“(O,x)  =  0,  0<ir<I 
and  inhomogeneous  boundary  conditions 

(2.6)  ^-*-(t,0+)  =  Ro<l>-(t,0-^),  </>-(!,£-)  =  Ri<i>+it,L-)  +  F(t) 

where  F  is  related  to  the  boundary  source  /  for  the  wave  equation  (2.4) 


(2.7) 


.^,((,0)  -  c(c^r(t,0)  =  0, 

1  +  (*,  <  >  0 


(2.8) 


t 

V)  =  s 


For  a  source  at  the  left  bormdary  (2.6),  (2.7)  are  to  be  replaced  accordingly.  The  first  order 
system  with  a  point  source  at  ^  e  (0,  L)  is 


a;)  +  (<,  x)  =  cF{t)6^, 

<f>~ (t,  x)  -  x)  =  cF{t)6^, 


0  <x  <  L,  0  <t  <T 


where,  as  in  (2.8),  2cF  is  the  primitive  of  the  point  source  /  in  the  wave  equation 


(2.10) 


<j>uit,x)  —  c^<l>xx{t,x)  =  /(t)^(,  0<x<L,  0<t<T. 


The  application  of  (dt  —  cdx)idt  +  cdx)  to  (j)'*'  +  ,  where  <j)'^ ,  <f>~  satisfy  (2.9),  (2.2),  (2.5) 

or  (2.1),  (2.6),  (2.5)  shows  that  these  systems  are  generalizations  of  (2.10),  (1.2)  or  (2.4), 
(2.7)  respectively,  with  zero  initial  conditions  =  V’o  =  0  in  (1.1.2). 

For  ^  €  (0,  L)  we  denote  by  the  space 

V(  =  {{•!>*, r)  6  H'ao.L)  \  («})'  1  ■^■^(0+)  =  J7<,r  (o+),r  (£■)  =  k,<a-(l-)} 

endowed  with  the  norm  on  (0,1-)  \  {0- 
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Deflnition  2.2.  Given  f  €  (0,i)  and  /  €  1^(0, T),  a  pair  (4,+  ,^-)  of  functions  4,*,  4,-  : 

[0  T)  X  (0,1)  E  wjlJi  4,+{0,x)  =  if-(0,x)  =  0,  0  <  I  <  L,  is  calied  a  soJutior 
a)  of  system  (2.1),  (2.6),  (2.5)  if 


(ii)  for  all  X  e  (0,L),  the  functions  [0,r]  R  de£ned  by  t  (f>-^(t,x),  t  h->  <l>-(t,x)  are 
absolutely  continuous  on  [0,T], 

(iii)  for  aJlte  (0,T),  (2.6)  holds, 

(iv)  for  aJi  t  G  (0,  T),  (2.1 )  holds  for  almost  all  x  E  (0,  L); 
b)  of  system  (2.9),  (2.2),  (2.5)  if 

(i)  the  mapping  [0,T]  V^,  defined  byt*-*  {(f>'^(t,-),4>~{t,-)),  is  in  C([0,T];  VJt), 

(ii)  for  all  X  E  (0,  L)  \  the  functions  [0,  T]  — »•  R  defined  byt*-^  x),  i 
are  absolutely  continuous  on  [0,  T], 

(iii)  for  all  t  E  (o,T).  .?i+(<,^+)  -  =  r(<,r)  - 

(iv)  for  all  t  E  (0,T),  (2.1)  holds  for  almost  all  x  E  (0,L)  \ 

Theorfim  2.2.  Let  f  E  L^(Q,T),  ^  E  (0,lt). 

a)  The  system  (2.1),  (2.6),  (2.5)  has  the  unique  solution 


<f>+(t,x)  =  f^Ro(RiRo)’'F  + 1)  x[L  -  ct,LK-2kL  -  x), 

k=o  \  c  / 

r  (*,  x)  =  f!„)‘F  (rSi±iM±£  +  ()  xir  -  cf,  L](-2kL  + 1). 

k=0  \  ^  / 

b)  The  unique  solution  of  (2.9),  (2.2),  (2.5)  is  given  by 


<^+(<,x)  = 


Y^{RoR^)^f( 

k=0 

oo  , 


^  —  2kL  —  x 


—  2kL  —  X 


+  x[C,C  +  ct]i2kL  +  x) 


2kL  —  x), 


(i>~{t,x)  =  '^{RiRo)’‘F  -  +  <")  x[^  -  ct,C]{-2kL  +  x) 

k=Q  \  ^  / 

+  f;  i2i  (i?oi2i  +  A  xlC,  e  +  ct]{2kL  -  x). 

k=i  \  ^  / 


Proof.  The  formulae  are  constructed  by  folding  the  solution  of  the  problem  on  all  of  R  (no 
boundaries)  into  the  domain  (0,L)  with  reflection  coefficients  Ro^Ri-  Note  that  for  any 
t  >  0  all  the  sums  sire  finite,  because  for  large  enough  fc  G  N  all  cheiracteristic  functions  are 
evaluated  outside  their  support.  Thus,  the  solution-properties  of  the  given  pairs 
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can  be  verified  by  considering  the  sums  term  by  term.  F  being  absolutely  continuous, 
the  regularity  with  respect  to  x  emd  t  follows  from  the  fact  that  F  is  zero  wherever  the 
characteristic  functions  have  time  dependent  steps.  That  the  norm  of  •)) 

varies  continuously  with  t  is  due  to  the  steady  propagation  with  finite  speed  c  (cf.  section  2.2). 
As  to  imiqueuess,  the  difference  of  two  solutions  of  one  of  the  inhomogeneous  problems 
(extended  continuously  at  ^  in  case  b))  is  a  solutioix  of  the  homogeneous  problem  (2.1), 
(2.2),  (2.5)  and  thus  equals  zero  by  Theorem  2.1.  □ 

For  (  €  (0,  £)  let  5^  denote  the  operator  that  maps  any  inhomogeneity  /  €  X^(0,T)  to 
the  restriction  to  (0,i)  \  {^}  of  the  spatial  derivatives  of  the  solution  of  (2.9),  (2.2),  (2.5), 
i.e. 

S(m = ■))  €  i,^((o,i)  \  {{))“. 

where  <f>~  are  given  in  Theorem  2.2  b).  The  point  ^  is  spared  out  to  exclude  ^-terms  in  the 
derivatives  of  •),  ^~(f,  Oi  whereas  the  velocity  potential  •)  +  •)  is  continuous 

at  Throughout  we  can  choose  S^f(t){^)  =  (0,0)  as  a  possible  extension  to  (0,X).  For 
notational  convenience  we  let  =  L,  =  0  and  denote  by  So  and  Sl  the  operators  that  map 
boimdary  sources  /  €  L^{0,T)  to  the  spatial  derivatives  •)>  ■))  solution  of 

(2.1),  (2.6),  (2.5),  where  (for  the  right  hand  boundciry)  are  given  in  Theorem  2.2  a). 

An  immediate  consequence  of  Theorem  2.2  is 

Corollary  2.2.  For  any  ^  €  [0,  L],  :  i^(0,  T)  -*  H  is  a  bounded  linear  operator. 

2.2.  Propagation  and  reflection.  The  components  of  the  solutions  given  in  section  2.1 
evolve  by  shifts  with  speed  c  in  positive  and  negative  direction  combined  with  inversion  and 
reduction  by  i2oi  -Ri  the  boundaries.  The  source  /  in  (2.10)  generates  a  wave  that  is 
symmetric  about  ^  €  (0,L).  More  specifically,  let  the  derivatives  of  the  two  components 
of  the  velocity  potential  (^‘•'(t,  •),  V’“(<,  •))  ==  be  given  at  some  time  <  >  0  and  let 

0  <  r  <  =  min(^,  L  —  0/^-  Then,  at  time  t  -{■  t 


(2.11) 


xl;+{t  +  T,x)  =  (^-7-  + 1  Xl^, (  +  cr](x) 

(  ll)'^{t,X  —  ct).,  CT  <  X  <  L 

\  —Roxl)~{t,CT  —  x),  0  <  X  <  CT 

t/»"(<  +  ^  /  (“T^  +t  +  r^x[(-  cr,^](x) 

J  V^“(t,x  -I-  cr),  0  <  X 

\  -Ritp'^{t,2L  -  X  —  ct),  L  -  c 


0  <  X  <  L  —  CT 

L  —  CT  <  X  <  L. 


This  is  to  be  understood  almost  everywhere  in  (0,  Z)  \  For  a  boundary  source  at  ^  = 
0(,f  =  L)  the  first  term  in  ip~{t  +  r,  x)  +  r,x))  is  to  be  deleted  with  —  tl  =  LJc. 

For  the  pair  ■=  S{(f>Q  ,(l)Q){t),  (2.11)  applies  for  0  <  r  <  Ljc  with  /  =  0. 

All  this  can  be  checked  by  differentiating  and  regrouping  the  formulae  of  section  2.1  and 
strict  adherence  to  the  convention  at  Theorem  2.1. 

2.3.  The  minimization  problem.  For  superposition  of  solutions  to  initial  conditions 
(^0  1  and  multiple  offending  sources  fi,-..  ,fn^^  points  , . . . ,  61  G  (0,  T]  define  g  ^  H 
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by 

n 

9(0  =  ,  ■#0  )(0  +  Y, 

i=l 

For  the  controllable  sources  located  at  r;i , . . . ,  G  [0,  L)  we  write 

u  =  (uu...,Um)eU:=  L\0,T)^ 
and  denote  by  B  the  bounded  linear  operator  U  —*  H 

m 

Bu{t)  =  '£s„Mt)- 

»=1 

In  this  notation  the  spatial  derivative  of  the  solution  of  the  first  order  system  that  corre¬ 
sponds  to  (1.1)  is  -f  Bu.  Rewriting  the  functional  (1.3),  note  that  for  solutions  {<p'^,<f>~) 
of  (2.1) 


for  all  t  €  (0,T)  and  almost  all  a:  €  (0,  L).  Therefore,  we  are  considering  the  following 
minimization  problem:  Given  {<f>o,<f>Q)  G  V  and  /i,...,/n  G  L^(0,T),  determine  u  £  U 
such  that 

(2.12)  J(ti)  =  min{  J(u) :  u  6  Cf}, 
where 

(2.13)  J{u)  =  {g  +  Bu,  Q{g  -|-  Bu))h  +  (u,  Ru)u. 

Here  Q  and  R  are  bounded  selfadjoint  operators  on  H  and  U  resp.,  defined  by 

=  />(g‘^(<,a:)<^+(<,x),g"(<,x)<^“(t,a:)) 
jR(ui,...,  u,„)(t)  =  (ri(<)ui  (<),...,  r,„(f)um(0)- 

The  two  components  of  the  waves  may  be  weighted  separately,  but  the  functionals  in  (1.3) 
and  (2.13)  coincide  if  =  q~  =  q.  We  make  the  assumption  that  the  weights  q"^  ,q~  G 
L'^(0,  r;i®®(0,L))  and  r,-  €  L®®(0,  T)  are  such  that 

(2.14)  R-fB*QR>0, 

where  B*  :  H  —*  U  is  the  adjoint  operator  of  B.  Then  the  unique  minimizer  of  the  quadratic 
functional  J  is  given  by  (eg.  [B],  (5.2.4)) 

u  =  -{R  +  B*QB)-^B*Qg, 
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3.  Discretization 


In  order  to  obtain  a  discrete  version  of  the  minimization  problem,  we  consider  a  grid  of 
charactenstics  of  (2.1)  in  the  (<,a:)-plane  with  knots  at  the  points  (tk,Xj),  where  Xj  =  jAx, 
j  =  0, . . . ,  iV  and  tk  =  kAt,  k  =  0,. . . ,  K.  We  assume  L  =  NAx,  T  =  KAt  and  Ax  =  cAt. 
Furthennore  all  sources  should  be  located  at  meshpoints,  i.e.  •  i^m}  C 

(xq  ,  .  .  .  ,  X ^  ) , 


3.1.  Piecewise  constant  sources  and  spline  initial  conditions.  Define  the  finite¬ 
dimensional  spaces 


Uk  =  {/  :  [0, T]  R  I  /  is  constant  on  [tjk-ijtjt),  fc  =  1, . . .  ,K], 

Vn  =  €  V  I  (f)'^ axe  first  order  splines  with  respect  to  {xj}^o}> 

Hn  =  ^  H  1 are  constant  on  [xj_i,Xj),j  = 

A  first  order  spline  with  respect  to  {xj}j^Q  is  a  continuous  function  [0,  i]  — ^  R  that  is  a 
polynomial  of  degree  one  on  [xj_i ,  x^],  j  —  1, . . . ,  JV.  Being  subspaces  of  X^(0,  T)  and  H,Uk 
and  Hn  axe  endowed  with  the  corresponding  L^-products.  For  the  A:th  basis  element  in  Uk 
we  choose  ^  —  1»  •  •  •  and  we  use  the  pairs  (x[xj_i  ,Xjf), 0)  and  (0, ^^j)) 

for  the  jih  and  {j  +  N)ih.  basis  elements  in  Hn,  ;  =  1, . . . ,  AT. 

Lemma  3.1.  If  0^)  €  V}v  and  fi  €  Uk,  i  =  1,. .  ■  ,n  then  gitk)  €  Hn,  A:  =  0, . . . ,  A . 
If  Ui  €  Uk,  i  =  I,. . .  ,m,  then  Bu{tk)  €  Hn,  k  =  0,...,K. 

Proof.  Because  of  Ax  =  cAt,  the  statements  of  the  lemma  follow  from  (2.11).  □ 

For  initial  conditions  €  Vn  and  sources  /,•  G  Uk  th®  coordinate  vector  p  of 

the  {K  +  l)-tupel  {g(to), . .  •  ,y(<fc))  e  has  2N{K  -b  1)  real  entries.  The  coordinate 

vector  u  of  u  =  (ui, . .  .,Um)  €  U'ff  is  in  R*”^.  The  operator  that  maps  any  u  G  U^  to 
(Sw(<o),  •  •  • ,  Bu{tK))  has  a  2N{K  -b  1)  x  mK  matrix  representation  B.  With  this  notation, 
the  time-sampled  spatial  derivative  of  the  solution  to  the  first  order  system  that  corresponds 
to  (1.1)  is  represented  by  the  vector  g  +  B_u. 


3.2.  Discrete  minimization.  We  replace  the  minimization  problem  (2.12)  by  the 
following  one:  given  (<^^ ,  <^^)  €  V}v  and  /i , .  •  • ,  /«  €  Uk,  determine  u  G  U^^  such  that 

(3.1)  J(u)  =  min{J(M)  :  u  G  }, 


where  J{u)  again  is  defined  by  (2.13). 
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With  the  numbers 


tk  j.. 

>j  ~  ^  J  J  q'^{t,x)dxdt,  =  p  J  j  q'^{t,x)dxdt, 


Pt 


tk-l  Xj-l 

th  Xj-c{t-tk-l) 


tk~l  X^-_i+c(t  — 
tk  XI 


Pk,j  =  P  J  J  q  {t,x)dxdt,  qk,j=P  J  j  ^  {t,x)dxdt, 


tk-l  *J-1 

tk 


tk-\  Xj-c(t-tfc-i) 


'i,k  =  J  ri{t)dt. 

tk-i 


we  define  the  2N  x  2N  matrices 


Qo  —  diag(5j^  J )  •  •  •  J  9i^7v  >  Pi  ,i » •  •  • » Pi  ,n)» 

Qk  =  diag(p]J‘  j  +  1, .  •  •  +  ^M-i.yv'Pifc+i,!  "I"  9jb,iv 

Qk  =  diag(pj^  J, . . .  ,ptc,Ni*lK,i^"  ■  ^ 9k, n)' 
and  the  square  matrices 


'Pjfc+1,N  +  9ik,jv)’ 
k  = 


g  =  diag(Qo,...,(5f('), 

R  =  diag(ri,i n, ,  r^.i ,  • . . ,  rm,K)- 


of  size  2N{K  +  1)  resp.  mK. 

A  minimization  problem  written  in  terms  of  these  matrices  is:  determine  u  G  such 
that 

(3.2)  Z(m)  =  inin{  J(u) :  u  6  R”*^}, 

where 

Hm)  =  (£  +  B_u)^Q{g  +  B  u)  +  u^Ru. 

We  assume  R  +  BT  QM.  >  0- 

Theorem  3.1.  The  coordinate  vector  u  of  the  solution  u  of  (3.1)  is  the  solution  of  (3.2). 

Proof.  Let  u  G  and  u  =  col(wi,i , . . . ,  ui^k-,  •  •  •  ?  Wm,i  j  •  •  •  >  Um,K)  he  its  coordinate  vector. 
Then 

m  K 

{u,Ru)u  = 

i=l  ifc=l 

Given  spline  initial  conditions  and  piecewise  constant  sources,  let  h{t)  =  {h'^(t,  •),  h~(t,  •))  = 
{g  +  Bu)(t)  and  denote  the  coordinate  vector  of  (h{to), ... ,  h(tK))  G  by  h  =  col(h^j, 


tk 

J  ri(t)u^i^dt  =  y^Ru. 
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•  •  • » ^t,N^  ^o,n  •  •  • )  ^o,Ni  •  •  • )  •  •  •  >  \i  -  •'  i^K  n)-  movement  of  the  two  com¬ 

ponents  of  h{i)  over  (xj_i,a;j)  during  tk-\  <i  <tk  is  indicated  in  Fig.  3.1. 


Xj  -c{t-ik-i) 


*j-i  +  c{t  -  tit -l) 


‘fc.i 


Xj-l 


Xi 


Fig.  3.1;  Piecewise  constant  components  at  time  t  6 


We  see  that 


tk  Xj 


j  j  +q~(t,x)h  {t,x)‘‘)dxdt 

= phi'Vi + 


tk~l  *i-l 


for  all  it  =  1, . . . ,  if  and  J  =  1, . . . ,  iV.  Using  this,  we  get 


K  N 


{KQh)H  =  ml  (g+(f,  x)h+{t,  xY  -t-  q-{t,  x)h-{t,  xf)dxdt 


^ 
N 


N 


=  S  +  PlJ^OJ  +  +  ^K,j^K,j 

j=l  J=1 

*:=!  i=l  k=2j=l 

=  h^Qh. 


Therefore,  J(u)  =  J(u).  □ 

Thus,  the  vector  representation  u  of  the  solution  to  (3.1)  solves  the  linear  system 
(3,3)  {E  +  R^QR)m.  -  -R^Qg- 
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4.  Numerics 


WiUiout  reference  to  the  particular  structure  of  system  (3.3),  its  solution  is  a  standard 
task.  However,  since  the  discretizations  in  time  and  space  are  coupled.  At  =  Ax/c,  the  size 
of  system  (3.3)  is  large  if  high  spatial  resolution  is  required  and/or  if  T  is  large.  Then  the 
number  of  entries  in  the  full  matrix  B  is  beyond  the  memory  capacity  of  standard  hardware. 
To  overcome  this  problem,  we  present  efficient  algorithms  for  the  computation  of  the  right 
hand  side  of  (3.3)  and  the  nonzero  entries  of  the  sparse  matrix  B'^QB. 

The  vector  £  can  be  produced  by  discrete  simulation  of  the  propagation  and  reflection  of 
the  components  •  2iV)  is  given  by  the  initial  conditions.  Instead  of  shifting  the 

arrays  that  hold  ,  (f>j^  at  each  time  step,  it  is  more  efficient,  in  particular  when  only  a  few 
of  the  time  steps  have  to  be  assigned  (Example  4.3),  to  keep  the  arrays,  except  of  replacing 
the  entry  for  <?ii'(0'*')  by  —Rq  times  the  previous  entry  for  <^j(0'^)  and  the  entry  for  <f)~{L~) 
by  ~-Ri  times  the  previous  entry  for  <f)^{L~).  To  accomplish  this,  an  integer  variable  for  the 
current  index  of  ^^(0^)  is  initialized  to  1,  decreased  by  1  at  each  time  step,  or  reset  to  N  if 
it  is  1.  The  same  boundary  index  is  applicable  for  both  components,  if  the  array  for  <f>~  is 
arranged  in  reversed  order.  To  keep  track  of  the  indices  where  to  add  the  contributions  of 
the  offending  sources  /j,  we  have  to  initialize  and  update  integers  that  point  to  the  source 
locations  i  =  1, . . . ,  n.  The  assignment  of  the  entries  of  <f>~  to  g{2kN  + 1  :  2(Jfc  +  l)iV) 
gives  the  kth.  group  of  £,  k  =  For  the  subsequent  multiplication  with  is  it 

convenient  to  store  the  groups  in  reversed  order.  Q  being  a  diagonal  matrix,  Q  ^  is 
obtained  from  g  by  2N{K  4-1)  multiplications. 

For  the  computation  of  —B^Q^  g,  consider  the  structure  of  Using  2c^  as  unit  for 
the  source  amplitudes  fi(t),  u,(<),  we  omit  the  factors  ^  in  (2.11).  Then  the  kth  basis 
elements  of  Uk  produces,  at  time  tkj  a  <j>~  pulse  to  the  left  and  a  <f>'^  pulse 

— iCj+:)  to  the  right  of  the  source  location  (boundary  sources  produce  only  one  such 
pulse).  These  pulses  travel  one  space  index  per  time  step  and  are  reflected  at  the  boimdaries. 
Thus,  the  fcth  row  of  a  source  block  in  B^  steirts  with  k  2iV- groups  of  zeros  and  continues 
with  K  —  {k  —  1)  2iV-groups  each  containing  (at  most)  two  nonzero  propagating  entries. 
Initisilizing  these  entries  to  q^l  according  to  the  source  location,  and  updating  their  indices 
and  amplitudes  according  to  propagation  and  reflection,  the  fcth  entry  of  the  ith  source  in 
—B^Qg  is  accumulated  by  adding  the  products  of  the  pulse  amplitudes  with  those  entries 
of  — Q  g  that  are  determined  by  the  updated  indices.  Repeating  the  entire  procedure  for 
i  =  1, . . . ,  m,  —  B^Q  g  is  computed  avoiding  the  zeros  in  (which  is  a  mK  x  2N{K  +  1) 
matrix  with  at  most  mK  x  2K)  nonzeros). 


Next,  we  describe  an  algorithm  to  compute  ^ QB_  for  time  independent  weights  q^{t^x) 
=  Then  Q  =  diag(lQi ,  Qi , . . . ,  Q, ,  iQi )  with  Qi  =  diag(g,'‘', 

•  •  •  1 9/v?  9r)  •  •  • » 9n)  ^  The  algorithm  is  beised  on  the  fact  that  two  waves  emerging 

from  two  sources  contribute  to  an  entry  in  B  B  only  if  the  pulses  that  are  generated  by 
the  basis  elements  of  Uk  overlap  during  their  travel  forth  and  back  in  the  duct.  Two  pulses 
overlap  starting  with  the  generation  of  the  later  one,  or  never.  Therefore  B  B  consists  of 
K  X  K  blocks  that  are  banded,  each  nonzero  diagonal  containing  the  accumulation  of  the 
products  of  two  overlapping  travelling  pulses. 


1.3 


An  entry  m  a  diagonal  can  be  computed  from  the  one  south  east  to  it  by  adding  the 
contribution  that  occurs  because  both  sources  are  active  one  cycle  earlier.  For  doing  this, 
It  IS  convenient  to  imagine  the  overlapping  pulses  tra\elling  backward  in  time,  adding  con¬ 
tributions  to  what  has  happened  later.  At  each  time  step,  the  product  of  the  overlapping 
pulses  has  to  be  multiplied  by  or  qj ,  according  to  the  current  direction  of  propagation 
and  the  current  location  of  the  overlapping  pulses  (by  half  the  value  of  qf  or  qJ  for  the  last 
time  step).  For  the  four  lower  diagonals  near  (resp.  on)  the  main  diagonal  of  the  block, 
the  location  of  the  entry  in  the  most  southern  row  is  determined  by  the  distance  of  the  two 
sources  of  the  block  and  their  distances  to  the  boundaries;  the  product  of  the  overlapping 
pulses  is  1,  i?o,  i?i,  RqRi  respectively.  These  four  diagonals  are  filled  up  one  after  the 
other:  The  elements  to  the  north  west  of  the  most  southern  entry  are  computed  until  the 
overlapping  pulses  arrive  at  the  boundary.  From  then  on  the  diagonal  entries  are  computed 
in  groups  of  N  elements  corresponding  to  the  propagation  back  and  forth  in  the  duct.  At 
the  beginning  of  each  group  the  product  of  the  pulses  is  multiplied  by  the  square  of  the 
appropriate  reflection  coefficient.  The  diagonal  is  complete  when  the  most  western  column 
of  the  block  is  reached. 

Each  time  a  nonzero  entry  to  one  of  these  four  diagonals  is  determined,  the  entries  within 
the  block  in  the  same  row  and  2j  N  columns  to  the  west  of  it  are  computed  by  multiplication 
with  {RoRiy.  These  represent  the  interaction  of  a  source  that  was  active  at  the  same 
location  but  j  full  periods  (=  i2iV  cycles)  earlier  than  the  one  just  considered  -  its  pulse 
was  reflected  2j  times  before  overlapping  with  the  later  pulse. 

Analogously,  the  diagonals  above  the  main  diagonal  of  the  block  are  built  up  from  south 
east  to  north  west  with  the  full  periods  north  of  the  current  entry. 

These  methods,  coded  in  Fortran  90  (NAGWare  compiler  1.2),  create  disk  files  that 
contain  the  right  hand  side  and  the  nonzeros  of  the  lower  triangle  of  the  coefficient  matrix 
in  (3.3).  The  latter  one  consists  of  lines  of  the  form  row  index,  column  index,  value  of 
the  entry.  The  files  are  loaded  into  Matlab  4.1  where  the  coefficient  matrix  is  set  up  as  a 
sparse  symmetric  matrix  A  and  the  right  hand  side  is  a  full  vector  rhs.  The  solution  u  is 
then  obtained  by  the  Matlab  command  u  =  A  \  rhs  that  invokes  sparse  matrix  arithmetics 
[M].  The  generation  of  the  optimal  waves  is  again  done  in  Fortran  90,  analogously  to  the 
generation  of  g.  Finally,  we  use  Matlab  for  the  graphics. 

In  the  following  examples,  t  is  given  in  seconds,  x  in  meters,  c  =  344.  The  relation  of 
the  diagonal  elements  of  Q,  denoted  by  q^{x),  to  the  weights  q^{t,x)  is  given  above.  We 
use  time  independent  functions  for  r,(<),  so  that  the  entries  of  R  eire  r,-  =  r,(<)A<.  The  step 
functions  w«(<)  are  measured  in  units  of  2c^.  Then  J(u)  =  J(u)  =  h^Qh  4-  u^R  u. 

Example  4.1.  L  =  10,  Ro  =  l.R^  =  0.5,  N  =  40,  K  =  137,  n  =  l,^i  =  4,m  =  1,t;i  = 
0,ri  =  0,q^{x)  =  1,  zero  initial  conditions.  This  introductory  example  is  comparable  to 
the  problems  in  [CNE].  The  offending  source  /i  =  x[^i97*2o)  emits  a  positive  and  a  negative 
pulse  at  time  <20  =  20Ax/c  ~  0.0145.  To  reduce  the  rightgoing  pulse  of  the  offending 
wave,  the  optimal  control  source  acausaly  emits  a  pulse  v/ave  at  time  <4:  ui{t)  =  —0.68  for 
i  £  [^3»^4)  •  To  annihilate  the  pulse  that  comes  from  rji  in  negative  direction,  ui{t)  =  —  1 
for  t  e  [<35, <30);  to  annihilate  the  reduced  pulse  that  is  reflected  at  the  right  boundary, 
u\{t)  =  —0.16  for  i  £  [<835 <84)-  After  that  <j)  is  constant,  silence.  The  3D  graphics  show  the 
velocity  potential  at  the  grid  points  {tk,Xj)  connected  by  straight  lines  that  are  parallel  to 
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Ui('  axes.  In  the  plane  of  the  t,x  axes  a  contour  plot  indicates  curves  of  constant  <i>{t,x). 
Compared  to  the  uncontrolled  process  J(0)  =  lie.!  2,  the  optimal  acoustic  (r,  =  0)  energy 
IS  J(h)  =  26.88. 

Example  4.2.  The  configuration  is  as  in  the  previous  example,  except  that  the  control 
source  now  is  located  in  the  interior,  t]\  =  2.  The  optimal  control  thus  generates  leftgoing 
pulses  so  that  the  annihilation  of  is  never  completed.  J{u)  =  56.95. 

Example  4.3.  L  -  IQ,Rq  =  Ri  =  0.2,  N  -  25,  K  -  148,  n  =  l,^i  =  10,  m  =  2,r?i  = 
3.2,  rj  =  0.5, 7^2  =  7.2,  r2  =  0.5,  g^(x)  =  x[0i3.2).  The  offending  source  at  the  right  bound¬ 
ary  runs  at  50  Hz,  fi{i)  =  sinlOOTrfjt  for  t  £  As  harmonic  initial  conditions  we 

take  (j)^ ,  <j>~  from  the  last  time  step  of  the  (uncontrolled)  generation  of  g  that  starts  form 
pre-initial  conditions  <^J(— 10,  •)  ~  <^“(—10,*)  =  0  with  f\  running  for  —10  <  t  <  0.  The 
optimal  control  causes  a  transient  (f  >  0)  toward  another  periodic  state  with  a  flat  velocity 
potential  on  x  (0,3.2).  Note  that,  because  of  r2  ^  0,  U2  switches  to  0  as  soon  as  its 

signals  cannot  influence  the  waves  in  (0,3.2)  any  more.  We  get  J(0)  =  659.46,  J{u)  =  20.65. 
The  flatness  of  (f)  in  (0, 3.2)  can  be  influenced  by  the  nagnitudes  of  ri ,  r2  relative  to  q^.  For 
our  choice  of  ri ,  r2  we  see  some  variation  of  <f>  throughout.  For  smaller  design  parameters 
ri,r2  the  oscillations  of  <f>  become  invisible  on  (<i8,<A')  x  (0,3.2).  However,  r2  =  0  would 
yield  a  singular  coefficient  matrix  since  q^{x)  =  0  on  (3.2, 10). 

Example  4.4.  To  test  the  applicability  of  our  methods  to  somewhat  more  realistic  acoustic 
events,  we  sample  the  first  2  seconds  of  a  wellknown  theme  of  classical  music,  /i(<)  = 
xlfi  ^^5^)sin3927rt  -f  x[li2)sin3127r<,  t  £  [0,2],  and  feed  it  to  a  source  at  =  2.5 
in  a  duct  of  length  L  —  5  with  reflection  coefficients  J2o  =  0,  i?i  =  0.7  and  zero  initial 
conditions.  We  seek  to  reduce  the  sound  near  the  boundaries:  q^{x)  =  x(0i  1)  +  x(4, 5).  Let 
two  control  sources  be  located  at  t/i  =  1,  r/2  =  4  weighted  by  n  =  r2  •=  0.  We  would  like  to 
have  spatial  resolution  of  about  30  grid  points  for  the  shortest  wave  length  of  the  offending 
wave,  344/196/30  =  0.059,  so  we  choose  Ax  =  0.05,  iV  =  100.  Then  T  =  2  =  KAt  leads  to 
K  =  2c/ Ax  =  13760. 

The  optimal  control  annihilates  the  components  of  the  wave  that  cross  t/i  ,  r/2  in  outgoing 
direction.  This  amounts  to  total  reflection  at  the  boundaries  of  (t/i,t/2),  which,  roughly 
speaking,  in  general  leads  to  increasing  sound  in  its  interior  where  the  offending  source 
is  located.  This  is  not  the  case  for  the  symmetric  configuration  chosen  here.  We  give 
graphics  for  the  transients  at  <  =  |  and  f  =  |.  At  the  end  of  the  first  tone,  the  graphics 
show  a  transient  toward  a  smoother  velocity  potential.  At  the  end  of  the  third  tone,  the 
controls  generate  waves  that  are  absorbed  by  the  offending  source,  so  that  ^  =  0  (although 
ri  =  r2  =  0  and  q^  =  0  in  (t/i,t/2))  until  the  fourth  tone  of  lower  frequency  begins  at  <  =  1 
(the  contour  curves  reveal  computational  inaccuracies).  Because  of  the  totally  reflecting 
behaviour  of  the  control  sources,  there  is  no  deimping  in  [t/i,t/2]  and  the  process  does  not 
converge  to  a  periodic  state  as  <  — >  2.  In  fact,  the  choice  of  ilo,  does  not  affect  the  optimal 
control  for  this  configuration;  however,  i^o  =  0  drastically  reduces  the  number  of  nonzeros 
in  (3.3).  The  full  matrix  B  would  require  over  600  gigabytes  of  storage  (8  bytes  per  entry), 
whereas  the  sparse  matrix  A  of  density  0.0003  has  2.4  megabytes  (including  the  memory  for 
the  indices).  The  sparse  arithmetics  to  solve  Au  =  rhs  require  less  tliEin  10.3  megaflops. 


15 


1 


5.  Approximation  and  Convergence 

We  again  assuire  that  L,  T  and  the  location  of  the  sources  are  such  that  there  exist  integers 
iV,  AT  €  N  with  •  ,r/m}  C  {jL/NJ  =  1,. . .  ,N  -  1}  and  K  =  cTN/L. 

For  ^  G  N  and  ^  L)  let  5/^  be  the  first  order  spline  interpolation  of  <j)  with  respect 

to  the  mesh  i.e. 


st(f>{3Lim)  -  <f>ULfeN\  j  =  o,...,eN. 

There  exists  a  co  >  0  such  that  for  all  ^  6  N  and  <!>  G  ^^(0,  L) 

(5.1)  ||s/<^  -  ^i|//i(o,L)  <  co||<^''||/^ 

(e.g.  [Sch],  Th.  2.5).  Moreover,  let  tt/  denote  the  orthogonal  projection  of  T^(0,T)  onto 
UiK^  i.e. 

kTjlK 

7r//(<)  =  ^  j  fdr  iorte[{k-l)T/eK,kT/eK),k  =  l,...,eK. 

(k-l)T/tK 

For  aU  /  G  L^{0,T)  (cf.  [BB],  p.  176) 

(5.2)  \\irtf  —  /l|L»(o,r)  -+  0  as  ^  >  co. 

Given  €  L^(0,T)  and  sufficiently  smooth  solution  u  of  (2.12)  is 

approximated  by  the  functions  u/,  that  are  the  solutions  of  (3.1)  wherein  )  and 

•  I'^tfn  axe  used  for  initied  conditions  and  sources  and  the  minimum  is  taken  over 
Ug^.  More  precisely  we  have 

Theorem  5.1.  If  (^^”,^0  )  G  V  fl  11^(0, 2)^  or  if  (<^^,<^^)  6  for  some  fo  €  N  then 
||ii  —  ti/lli;  — >  0  as  f  -♦  00. 

n 

Proof.  With  gt  =  S(si<^Q,si<f>o)  +  X)  the  optimal  control  u/  G  is  characterized 

»=i 

by  the  variational  equation 

a{{i(,u)  =  bi{u)  for  ail  u  G 

where  a  :  U  x  U  —*  R,  bg  :  U  R  are  symmetric  bilinear,  resp.  linear,  continuous  forms 
given  by 

a(u,  v)  =  (u,  {R  +  B*QB)u)u 
btiu)  -  -{B*Qgt,u)u. 

This  is  so,  because  u/  minimizes  the  functional  Jt{u)  =  a{u,  u)  —  2be{u)  +  {gt,  Qgi)H  in 
(cf.  [C],  Theorem  1.2;  the  constant  term  is  included  here  for  consistency  with  (2.13),  but  it 
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doe.',  e,  .  affect  the  minimization  problem).  Similarly,  u  minimizes  J{u)  =  a(u,u)  —  2b{u)  -f 
{fliQg'lH  in  U  with 

biu)  =  -{B*Qg,u)u 

and  thus  it  is  characterized  by 


a{u,u)  =  6(u)  for  all  u  G  f/. 

By  assumption  (2.14)  a  is  CZ-elliptic  and  we  can  apply  the  first  Strang  lemma  ([C],  Th.  26.1); 
there  exist  constants  c\ ,  C2  that  are  independent  of  such  that 

llw  -  ui\\u  <  Cl  inf{||u  -  u\\u  :  u  e  U^}  +  C2  sup{|5(u)  -  6/(u)|/||u||(;  ;  u  €  U^] 

The  first  term  on  the  right  hand  side  tends  to  zero  by  (5.2).  The  second  term  is  estimated 
using  Schwarz’  inequality  and  the  boundedness  of  B*Q 

\{B*Q{g-gt%u)u\/Mv  <  \\B*Q\\\\g  -  giU 


where 


n 

Hi?  -  wll//  <  ||•S'||||(^o’,^^)  -  ^  ||•S'^,•||||/<  -  ff</i||L»(0,T) 

i=l 

so  that  (5.1),  (5.2)  imply  the  result.  In  case  (^J,  4>o)  ^  ^toN  we  have  si<f>Q  =  <f>Q,  st<l)Q  = 
for  f  >  £o.  □ 

The  second  case  in  Theorem  5.1  is  applicable  when  the  problem  is  broken  up  into  pieces 
over  subintervals  of  [0,T].  But  this  is  not  equivalent  to  minimization  over  [0,T]  at  once. 

How  do  the  discrete  controls  u/  perform  when  applied  to  the  original  data  {<f>Q,<f>o)  € 
V  n  H^{0,Ly  and  sources  /,  G  L^(0,  T)?  The  answer  (in  the  sense  of  convergence)  is 
contained  in  the  theorem:  J{ut)  — »  <^(w),  since  J  :  U  — »  R  is  continuous. 

5.  Conclusions 

By  decomposition  of  the  one-dimensional  wave  equation  with  point  sources  zind  pointwise 
reflecting  boundary  conditions  and  by  appropriate  discretization  in  time  and  space,  the 
minimization  problem  with  piecewi.se  constant  sources  was  reduced  to  the  solution  of  a 
sparse  linear  system.  It  was  shown  that  the  waves  that  are  generated  by  discrete  forcing 
functions  and  initial  conditions  are  sufficiently  regular  so  that  the  problem  is  wellposed.  Four 
examples  demonstrated  the  applicability  of  algorithms  that  efficiently  set  up  the  sparse  linear 
system;  such  examples  provide  insights  into  the  optimal  control  of  waves  in  ducts.  Finally, 
it  was  proved  that  the  solutions  of  finite  dimensional  discretizations  converge  to  the  solution 
of  the  minimization  problem  with  square  integrable  sources. 
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